"""
生成命令
"""
import re
import os
import pandas as pd

commands = []

bed_df = pd.read_csv("./A188.Insertion.test.bed", sep="\t", header=None)

for index, row in bed_df.iterrows():
    item_commands = []
    chr_id = row[0]
    start = row[1]
    end = row[2]

    sv_dir = f"./{chr_id}_{start}"
    os.makedirs(sv_dir, exist_ok=True)

    # step1: 提取数据
    item_commands.append(f"""
        bcftools view --regions {chr_id}:{start}-{end} \
        /public/home/yunlzhang/data/Jilin/3613/reference_panel/JL.3613.imputed.vcf.gz \
        -Oz -o {sv_dir}/SV.vcf.gz
    """)
    # step2: 格式化数据
    item_commands.append(f"""
        zcat {sv_dir}/SV.vcf.gz | awk 'BEGIN{{OFS="\\t"}} /^#/ {{print $0; next}} {{$3=$1"_"$2; print $0}}' | bgzip -f > {sv_dir}/SV.id.vcf.gz
    """)
    # step3: plink 转换数据
    item_commands.append(f"""
        plink --vcf {sv_dir}/SV.id.vcf.gz --allow-extra-chr --geno 0.5 --maf 0.05 --recode --out {sv_dir}/SV.plink
    """)

    # step4: 生成haploview入口文件
    item_commands.append(f"""
        sed -i 's/-9/0/g' {sv_dir}/SV.plink.ped && awk 'BEGIN {{OFS="\\t"}} {{print $2, $4}}' {sv_dir}/SV.plink.map > {sv_dir}/SV.plink.info
    """)

    # step5: 执行haploview
    item_commands.append(f"""
        java -jar ~/software/Haploview.jar \
            -nogui \
            -pedfile {sv_dir}/SV.plink.ped \
            -info {sv_dir}/SV.plink.info \
            -minMAF 0.05 \
            -missingCutoff 0.5 \
            -maxDistance 5000 \
            -blockoutput \
            -blockCutHighCI 0.8 \
            -pairwiseTagging \
            -svg \
            -out {sv_dir}/haploview
    """)



    format_commands = "&&".join([re.sub(r"\s+", " ", command) for command in item_commands])
    commands.append(format_commands)

with open("commands.sh", "w") as f:
    f.write("\n".join(commands))
